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INTRODUCTION 

In semitransparent materials where thermal radiation can 
affect internal temperature distributions, transient behavior 
has been studied much less than steady-state. To obtain 
transient solutions, numerical procedures such as finite 
difference and finite element methods have been used to solve 
the radiative transfer relations coupled with the transient 
energy equation. Some of the literature has been reviewed in 
[1-3]. In [3] transient solutions were obtained for a layer 
with a refractive index larger than one with external con- 
vection and radiation at each boundary; these results using 
the exact equations of radiative transfer will be used for 
comparison with the present two-flux calculations. 

Various multi-flux methods have been discussed [4] as a 
simplification for computing the radiative flux term in the 
energy equation. For the general boundary conditions of 
external convection and radiation on a layer with diffuse 
interfaces it was shown in [5] that the two-flux method can 
be used to predict accurate steady-state temperature dis- 
tributions and heat fluxes. The purpose of this note is to 
show that the two-flux method can be used to obtain tran- 
sient solutions in materials with large refractive indices that 
are typical of ceramics. Good predictions of transient tem- 
perature distributions are obtained as verified bv comparison 
with implicit numerical solutions from [3] where the exact 
radiative transfer equations were used. The information on 
transient behavior is ol interest for computing heat transfer 
performance in ceramics being developed for high tem- 
perature use in advanced aircraft engines. 

An advantage of the two-flux method is that isotropic 
scattering is included without any additional complication. 
Some transient results with large scattering are given to illus- 
trate scattering effects ; the solutions in [3] are for absorption 
only. 


ANALYSIS 

Energy and two-flux equations 

A plane layer of thickness D . Fig. 1, is a heat conducting, 
gray emitting, absorbing, and isotropically scattering medium 
with n ^ I, and its boundaries are assumed diffuse. The layer 
is initially at uniform temperature 1\ and is placed in sur- 
roundings so each boundary receives radiative energy and is 
subject to convection. Transient temperature distributions 
are to be obtained in the layer until steady-state is reached 
corresponding to the external radiation and convection con- 
ditions. 

The transient energy equation in dimensionless form is [3] 


Properties are assumed independent of temperature. The 
gradient of the radiative flux. cq T (X, r)/T.V, is obtained from 
the two-flux relation using the Milne Eddington approxi- 
mation [4, 6], 

cq Y ( .V, r) 

= * 1 ) ( 1 - Q ) [4/r t 4 ( A', T) - G{ a; t)] (2) 

where (7(>\,r) is related to <y r (A, t) by the equation, 

c(7(A, t) 

f~X' ~ ~ r). (3) 

The q f and (i are related to the positive and negative radiative 
fluxes shown in Fig. 1 by fliX.x) = q, [ (A,r )-q r (A'.t) and 
G(X,t) - 2 [ij; ( X, t ) -H/ r ( A , r)]. 


Boundary and initial conditions 

The convective boundary conditions on the sides of the layer 
are 


?t 

M v 


//, 

4A ,: 


[/,i ~/(0,D] 


(4a) 



Fig. 1. Geometry and nomenclature for transient radiation 
and conduction in a semitransparent layer with isotropic 
scattering. 
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NOMENCLATURE 


a absorption coefficient in layer [m ’] 

B the quantity 3h.f>( 1 -H) 
t specific heal of radiating medium 

[J kg 1 K ’] 

D thickness of semitransparent layer [m] 

O' the flux quantity 2 (#/ +q t ) [W m : ] ; 

(i — GioTf 

/j , , A , convective heat transfer coefficients at 

.v = 0 and D [W m K '] 

//,. //, dimensionless parameters. h,aTl and 

k thermal conductivity of layer [W m 1 K '] 

n refractive index of layer 

:V conduction radiation parameter, kj 4<t 7’, 3 D 

(f t radiative llux in the v-direclion [W m ] ; 

‘/l = ‘/l nT < . . 

V,i- i external radiation fluxes incident at a = 0 
and D [W m : ] 

ij ri , dimensionless radiation fluxes, 

qnlaTt.qJrf't 

T absolute temperature [K] ; / ~ T; 1\ 


gas temperatures for convection at a = 0 
‘ and D [K] 

f Fl , if dimensionless gas temperatures, / ^ I „ 

TfjT, 

T, initial uniform temperature of layer (used 
as a reference T) [K] 

v coordinate in direction across layer [m] ; 

X = v D. 


Greek symbols 
0 time [s] 

K n optical thickness of layer, {a + rxj D 

p density of semitransparent medium 

[kg m ■ ] 

p\ p° internal and external reflectivities at a 
boundary 

c t Stefan Boltzmann constant [W m ‘ K 4 ] 

<t s scattering coefficient in layer [m ! ] 

r dimensionless time, (4 oT\ fpcD)0 

Q scattering albedo, (j s , (a + <r s ). 


ft 

n 




■^1. 


(4b) 


The radiative boundary conditions must now be specified 
including the effects of internal and external reflections at 
the surfaces. By considering the incident and reflected fluxes 
at an interface the following boundary relations between G 
and were developed at each boundary [5] 


f7(0.x) — 4 1 — </,, - t7,(0,r) (5a) 

1 - p' 1 p l 

(7(1, x) =4-- — q {2 + 2 ]+f) </,(!. t). (5b) 

\-p l 1-// 

To begin the transient solution of equation ( l ) the specified 
initial condition is a uniform temperature 7 ’(a',0) = T { so 
t{X, ()) = 1. Initial distributions are also needed tor 0) 
and J(.V.O). By differentiation, equations (2) and (3) can be 
combined to eliminate either </ t (A,t) or (/(A.x) to give a 
second order equation tor either of these quantities. With 

t = l initially, these equations a re solved analytically to give. 

(7(A\0) = r,e'"' + C\e ' " v +4/r (6a) 

0) — (6b) 

AK,> 

where B = 3xf,(l -Q). C, and C_ are integration constants 
that are obtained by applying the boundary conditions (5a) 
and (5b) to equation (6a). The following quantities are 
defined : 


1 - p“ 

S> = — 4m" +4— - — i hi ■ 
1 -P ! 

Then the integration constants are. 

-pSj+tSi 

( ' -~-l*y+xS 
*S 2 -ySi 
C; ~~py+xS' 
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(8a) 


(8b) 


Numerical solution 

Starting with the initial /(A, 0) and y r (A,0) relations, equa- 
tion (1) was integrated forward in time using the following 
explicit finite-difference algorithm at the interior grid points : 


t(X, t + At) - f(X, t) + 


NAx 

(Axr 


[tiX+AX'i) 


- 2/(A\ t) + t(X—AX, t)1 


At vq f 

4 (TV 


(9) 


The r(x + Ax) at the boundaries were then evaluated using 
equations (4a) and (4b) with a three-point difference 
approximation for the temperature derivative 


r(0, t + At) 


11, AX 
~2N 


tf T 4/( AA\ r + At) - /(2AAT x + Ax) 


II, AX 

3+ ”2jv~ 


( 10a) 


a = 

1 + />' 

1-2 • 


(7a) 

r(l,x + Ax) 


1 -p' 

3k„ 


I LAX 


i if 

v* 

3h'i> 

(7b) 

_;+T: + 


AA', t + At) — /(I -2AA\x + Ax) 

H AX 

3+-frr- 


(10b) 


'/ = P e s B 

<) = K‘ ' H 

1 - p" 

, = -4/r+4 (/ r i 

1 -p l 


(7c) After advancing t{X) each time increment the radiant flux 
gradient in the last term in equation (9) must be advanced 
(7d) tQ t + At. This was done by solving equations (2) and (3) 
simultaneously for G(X, t + At) and r/ r (A\ x + Ax) using 
(7e) /(A, t + At) on the right-hand side of equation (2). Then 

<Ty r (A\ t + At)/TA' was evaluated from equation (2). The sim- 
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ultaneous solution of equations (2) and (3) was carried out 
using a fourth-order Runge-Kutta method with a shooting 
procedure to satisfy the boundary conditions at X = 0 and 
1. To begin the solution the value of q r ( X = 0) from the 
previous time step was used as an estimate, and the boundary 
condition equation (5a) was solved for G(X — 0). The solu- 
tion was then carried out by Runge-Kutta integration from 
X = 0 to X = 1 . The values of q t and G obtained at X = 1 


were checked to see if they satisfy the boundary condition in 
equation (5b). An iteration was performed on q r {X = 0) until 
equation (5b) w r as satisfied; the type of iterative method 
used is described in [7]. The shooting method used here is 
convenient for absorption optical thicknesses, aD < 8. This 
two-point boundary value solution method becomes difficult 
when there is a large aD that causes the conditions at the 
two boundaries to become less directly related. It is possible 





0. 0 0. 2 0. 4 0. 6 0. 8 1.0 

Dimensionless coordinate, X = x/D 

Fig. 2. Two-flux and exact results for transient temperature distributions in a layer initially at uniform 
temperature after exposure to radiation on one side and convective cooling on the other side ; no scattering, 
Q = 0. Parameters: N = 0.1, </ rl = 1 5 4 . <y r2 = 0.5 4 . //,=(). H 2 = 1. ^ = 0.5. (a) Optical thickness. 
= ().5 ; (b) optical thickness, k„ = 2 ; (c) optical thickness. K* n = 5. 
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that other numerical techniques could partially eliminate this 
difficulty: a method using a Green's function is presently 
being developed. Using the t{X) and G(X) the radiant flux 
gradient was evaluated from equation (2). The temperature 
distribution was then advanced to the next time increment 
using equations (9) and (10). 

After checking various grid sizes it was found that 41 
evenly spaced points across the layer gave accurate solutions. 
Corresponding to this grid size (A,V — 0.025) and for 
V = 0.1, as used for the results given here, the At for a stable 
explicit calculation was estimated from the criterion for solv- 
ing the transient heat conduction equation. The value 
Ar = 0.0025 provided stable solutions for all of the results 
calculated here. 


RESULTS AND DISCUSSION 

The transient temperature distributions given here start 
from a uniform initial temperature T(X. 0) = / so that 
/(TO) = 1. Figure 2 shows typical comparisons of two-flux 
results for /(.V. t) with those from [3] using an implicit finite- 
difference method, and the exact equation of transfer to 
evaluate the radiative flux gradient in the energy equation. 
The layer is heated on the hot side fT = 0) by a radiative flux 
equal to that from blaekbody surroundings at T sl = 1.57',; at 
\ = I, 7T = 0.5 T x so there is a net radiative cooling at that 
side. These illustrative results examine the thermal behavior 
of a layer that is convectively cooled only on the side away 
from where the radiative heating occurs. This simulates pos- 
sible conditions for the wall of a combustion chamber w here 
there is radiative healing from combustion gases on one side, 
and that side is not being film cooled. 

The three parts of Fig. 2 arc for optical thicknesses of 0.5. 
2 and 5 w ith no scattering. For K n — 0.5 the layer is somewhat 
optically thin. For K n = 2 the optical thickness is such that 
maximum internal radiative effects are expected ; for k,, — 5 
the layer is somewhat optically thick. Each part of the figure 
shows results for n = 1 and 2. The two-flux curves are solid 
or long dashes ; the numerical results using the exact transfer 
equations are medium or short dashes. When t — 1.5 the 
temperatures are within 1% of steady slate. At X — 0 the 
temperature profiles have a zero derivative from the absence 
of convective cooling at that boundary. The convective cool- 
ing at X = I produces a rapid temperature decrease near that 
boundary. 

The results using the two-flux method agree within reason- 
able error with predictions using the exact transfer equations. 
The largest deviations, which are for = 0.5. are only a few 
per cent, and agreement is much better for k d — 5. As // 
increases, internal reflections make the temperature dis- 
tributions more uniform. In most instances agreement of the 
two-flux results was a little better for n = 1 than for n = 2. 

The effect of scattering is illustrated in Fig. 3 for n = 1 and 
2. The optical thickness is constant. K n = 5. so an increase 
in scattering corresponds to a decrease in absorption. The 
result is that the transient temperatures are decreased with 
increasing Q. For // = 2 in Fig. 3(b) the temperatures are 
somewhat more uniform than for n = 1 in Fig. 3(a) [note 
that the ordinate scales are different in Figs. 3(a) and 3(b)]. 
C ompared with Fig. 3(a), increasing Q in Fig. 3(b) does not 
have as large an effect in reducing the temperatures. For 
n — 2 the layer has internal reflections that make scattering 
more effective in augmenting absorption. For SI = 0.99 this 
makes the temperatures larger for ;; = 2 than for n = 1 

CONCLUSIONS 

The two-flux method was used to obtain transient solu- 
tions for a plane layer including internal reflections and scat- 
tering. The layer was initially at uniform temperature, and 
was heated or cooled by external radiation and convection 
The two-flux equations were examined as a means for eva- 


luating the radiative flux gradient in the transient energy 
equation. Comparisons of transient temperature dis- 
tributions using the two-flux method were made with results 
where the radiative flux gradient was evaluated from the 
exact radiative transfer equations. Good agreement was 
obtained for optical thicknesses from 0.5 to 5 and for refrac- 
tive indices of 1 and 2. Illustrative results obtained with the 
two-flux method demonstrate the effect of isotropic scat- 




Fig. 3. Effect of scattering on temperature distributions in 
a layer initially at uniform temperature after exposure to 
radiation on one side and convective cooling on the other 
side. Parameters: k d = 5, N = 0.1, <y r , = 1.5 4 , ij, 2 = 0.5\ 
//, = 0. H : = F t g2 = 0.5. (a) Refractive index, n = 1 ; (b) 
refractive index, n — 2. 
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tering coupled with changing the refractive index. For small 
absorption with large scattering the maximum layer tem- 
perature is increased when the refractive index is increased. 
For larger absorption the effect is opposite, and the 
maximum temperature decreases with increased refractive 
index. 
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